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Abstract 

A common assumption in satellite imager-based cirrus retrieval algorithms is that the 
radiative properties of a cirrus cloud may be represented by those associated with a specific 
ice crystal shape (or habit) and a single particle size distribution. However, observations of 
cirrus clouds have shown that the shapes and sizes of ice crystals may vary substantially 
with height within the clouds. In this study we investigate the sensitivity of the top-of- 
atmosphere bidirectional reflectances at two MODIS bands centered at 0.65 pm and 2.1 1 
pm to the cirrus models assumed to be either a single homogeneous layer or three distinct, 
but contiguous, layers. First, we define the single- and three-layer cirrus cloud models with 
respect to ice crystal habit and size distribution on the basis of in situ replicator data 
acquired during the First ISCCP Regional Experiment (FIRE-II), held in Kansas during the 
fall of 1991. Subsequently, fundamental light scattering and radiative transfer theory is 
employed to determine the single scattering and the bulk radiative properties of the cirrus 
cloud. Regarding the radiative transfer computations, we present a discrete form of the 
adding/doubling principle by introducing a direct transmission function, which is 
computationally straightforward and efficient, an improvement over previous methods. 

For the 0.65 pm band, at which absorption by ice is negligible, there is little 
difference between the bidirectional reflectances calculated for the one- and three-layer 
cirrus models, suggesting that the vertical inhomogeneity effect is relatively unimportant. At 
the 2.1 1 pm band, the bidirectional reflectances computed for both optically thin (x = 1) 
and thick (x = 10) cirrus clouds show significant differences between the results for the 
one- and three-layer models. The reflectances computed for the three-layer cirrus model are 
substantially larger than those computed for the single-layer cirrus. Finally, we find that 
cloud reflectance is very sensitive to the optical properties of the small crystals that 
predominate in the top layer of the three-layer cirrus model. It is critical to define the most 
realistic geometric shape for the small “quasi-spherical” ice crystals in the top layer for 
obtaining reliable single-scattering parameters and bulk radiative properties of cirrus. 
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1. Introduction 


Cirrus clouds located in the upper troposphere and lower stratosphere are important 
to the Earth’s climate [Liou, 1986; Stephens et al., 1990]. They reflect solar radiation, 
absorb the thermal emission from ground and lower atmosphere, and reemit infrared 
radiation to space. The Moderate Resolution Imaging Spectroradiometer (MODIS) [King et 
al., 1992] on the recently launched Terra spacecraft will enhance the capability for 
monitoring cirrus clouds in comparison with previous generations of satellite instruments 
such as the Advanced Very High Resolution Radiometer (AVHRR). MODIS has a total of 
36 bands for studies of land, ocean, and the atmosphere, including the 1 .38 pm band for 
cirrus detection and correction [Gao and Kaufman, 199 5] , 0.65, 1.66, and 2.11 pm bands 
for implementing cloud microphysical and optical property retrieval techniques [King et al., 
1997], and the 8.5, 11, and 12 pm for applying infrared cloud property retrieval techniques 
[Ackerman et al., 1990]. MODIS data will allow improved retrievals of cirrus optical and 
microphysical parameters such as optical thickness of cirrus and mean effective size of the 
ice crystals in these clouds. 

Various algorithms have been developed to retrieve cirrus optical and microphysical 
properties in the past twenty years [e.g., Liou et al., 1990; Ou at al., 1993; Minnis et al., 
1993a,b; Rossow and Lacis, 1990]. They can be categorized into the techniques based on 
either infrared emission or solar reflection. The representative algorithm of the former is the 
method developed by Inoue [1985] for determining the infrared emissivity of cirrus clouds 
on the basis of the brightness temperature difference between 11 and 12 pm wavelengths 
with an assumption of implicit mean particle size. An infrared trispectral algorithm using 
8.52, 1 1, and 12 pm bands [Ackerman et al., 1990, 1998; Strabala et al. 1994] with a recent 
improvement [Baum et al. 2000a,b] form the basis of an infrared retrieval algorithm using 
MODIS infrared channels. The representative retrieval algorithm based on solar reflection is 
that developed by Nakajima and King [ 1990], who used a trispectral (0.75, 1.6, and 2.2 pm) 
method to simultaneously retrieve the optical thickness and mean effective particle size for 
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water clouds. This approach has been applied to the retrieval of the optical thickness and 
mean particle size of ice crystals for cirrus [Wielicki et al, 1990]. 

To develop a reliable retrieval algorithm for cirrus optical and microphysical 
properties, it is critical to generate reliable pre-calculated look-up tables of bidirectional 
reflectance for cirrus clouds over a practical range of effective sizes and optical thicknesses 
for various solar- and view-angular geometries. At present, most algorithms for retrieving 
cirrus optical thickness and effective size assume that all ice crystals in cirrus are of one 
specific habit, such as spheres, hexagonal plates, hexagonal columns, or fractal polycrystals. 
In addition, it is commonly assumed that a single size distribution is sufficient to determine 
the scattering properties of the ice crystals within the cirrus layer. However, the observations 
based on aircraft-borne two-dimensional optical cloud probe (2D-C) and balloon-borne 
replicator measurements [e.g., Heymsfield and Platt , 1984; Amott et al, 1994; Mitchell et 
al., 1996a,b; McFarquhar and Heymsfield, 1996, 1997] have demonstrated that a significant 
number of ice crystals in cirrus clouds are bullet rosettes, solid and hollow columns, plates, 
and irregularly shaped aggregates. In addition, Heymsfield and colleagues have showed that 
ice crystal habits and size distributions are vertically inhomogeneous in cirrus clouds. 
[http://box.mmm.ucar.edu/science/ciiTus/projects/FIRE]. 

Since satellite-based retrieval techniques essentially compare library computations of 
bidirectional reflectances to actual measurements in their implememtations, it is necessary to 
assess the effect of the inhomogeneity of the ice crystal sizes and shapes within cirrus on 
the radiative transfer calculations for generating the reflectance libraries. Our objective is to 
understand the effect of vertical inhomogeneity in the structure of cirrus clouds on their 
radiative properties. We employ fundamental scattering and radiative transfer theory to 
investigate the bidirectional reflectance of cirrus clouds at MODIS 0.65 and 2.11 |xm bands 
using in situ crystal habit and size distribution for two typical cases of midlatitude cirrus 
cloud systems. 
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Section 2 describes the data and models used in this study. In Section 3, we first 
define three-layer and single homogeneous cirrus models. We th$n derive the single- 
scattering properties associated with the two cirrus models. In section 4, we present the 
differences between bidirectional reflectances computed for the three-layer cirrus model and 
its one-layer counterpart. Also presented in this section is the sensitivity study regarding the 
shape effect of the small “quasi-spherically” nonspherical ice crystals on cloud reflectance. 
Conclusions are given in section 5. Finally, in Appendix we present a numerically stable 
radiative transfer model based on the adding/doubling principle expressed in a discrete form 
for calculating reflected and transmitted intensities resulting from multiple scattering and 
absorption of cirrus clouds. 

2. Data and Models 
a. Data 

The size distributions and ice crystal habit information obtained from two cases of 
replicator measurements carried out in Kansas during FIRE-II are used in this study. The 
balloon-borne ice crystal replicators were launched at 1337 UTC on 25 November 1991 and 
at 2045 UTC on 05 December 1991 and ascended at a rate of 4 m s' 1 while passing through 
the cloud layers. Cloud top temperature on 25 November was -57°C and on 5 December 
was -65°C, as measured by a radiosonde connected to the balloon package. The replicator 
collects particles in a liquid plastic solution, which is coated on a moving, 35--mm-wide 
transparent leader tape. The particles are preserved when the solvent in the solution 
evaporates, giving detailed ice crystal impressions and size spectra of crystals down to 10 
microns and below. The particles do not breakup upon impact on the replicator tape, 
because of the slow rate of ascent of the balloon. The efficiency with which the replicator 
collects small crystals has been quantified theoretically and experimentally [Miloshevich 
and Heymsfield, 1997], and concentrations used in this study are adjusted to account for the 
non-zero collection efficiencies of small particles. The ice crystals were sized throughout the 
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depth of the cloud, which, when their concentration was appropriately considered, led to 28 
size spectra in the vertical on November 25 and 33 spectra on December 5, each spectra 
representing approximately 100 m in the vertical. 

b. Single-Scattering Properties 

Ice crystal sizes in cirrus clouds fall predominately within the applicable size 
parameter regime of the geometric optics method at visible and near-infrared wavelengths. 
Each ray can be localized on the wavefront of incident radiation (electromagnetic wave); 
consequently, Snell’s Law and Fresnel’s formulas can be applied to trace the ray 
propagation and the electric field magnitude as well as polarization configuration associated 
with the ray. In the conventional geometric optics approach for deriving the scattering 
properties of a particle, the scattered field in the radiation zone is regarded as the 
superposition of diffracted rays and Fresnelian rays. The contribution due to ray diffraction 
can be solved using the standard Fraunhofer theory whereas the contribution of Fresnelian 
rays can be solved by a ray-tracing calculation. This approach suffers from several 
shortcomings, as noted by Yang and Liou [1995]. In particular, it is assumed that the 
extinction efficiency is 2 regardless of particle size. 

To overcome the shortcomings of the conventional geometric optics, Yang and Liou 
[1996] developed a geometric-optics-integral-equation approach (hereafter referred to as 
GOM2). According to fundamental electrodynamics theory, the scattered far field can be 
obtained if the tangential components of the electric and magnetic fields on the particle 
surface are specified. In principle, GOM2 employs the ray tracing technique to solve the 
near field on the particle surface and then map the near field to the far field via the following 
rigorous relationship: 

s 
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where n s and ? are unit vectors along the normal directions of particle faces and scattering 

direction, respectively, and k = In I A in which A is the wavelength of incident wave in 
vacuum. It should be pointed out that the partition of the diffraction and Fresnelian 
contributions to the far field is not required in GOM2. A simplified algorithm for GOM2 is 
employed to reduce the computational cost. 

To compute the extinction and scattering cross sections of ice crystals, we use the 
ray tracing technique coupled with the following exact electrodynamic relationships: 
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where E 0 is the electric field associated with incident wave, e= e r +/e, is the complex 


permittivity, and the asterisk indicates the complex conjugate operation. The volume 
integrals in Eqs.(2a) and (2b) are carried out along individual ray paths inside the particle 
via a Monte Carlo/ray-by-ray algorithm [Yang and Liou, 1997] that is a generalization of 
the well-known anomalous diffraction theory [van deHulst, 1957]. 

Because of the complicated mechanisms involved in ice crystal growth, such as 
sublimation, riming or aggregation of particles observed in anvils, the surfaces of ice 
crystals may be rough. Roughened ice crystal surfaces were observed in laboratory 
experiments [Cross, 1968] and also verified by recent in situ observations of tropic anvil 
cirrus clouds [Heymsfield et al., to be available at http://box.mmm.ucar.edu/science/cirrus ]. 
The scattering phase functions of roughened particles display less variation [Macke et al., 
1996; Yang and Liou, 1998] than their counterparts for ice crystals containing smooth 
facets. As articulated by Mishchenko et al. [1996] on the basis of ground-based 
nephelometer and aircraft radiance measurement of cirrus clouds [Foot, 1988; Francis, 
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1995; Gayet, 1995; Posse and von Hoyningen-Huene 1995], the scattering phase functions 
of cirrus clouds can be rather featureless with no appreciable halos. 

Based on these studies, we account for surface roughness specifically for ice crystal 
aggregates in this study. In the numerical computation, the particle surface is regarded as a 
number of small facets whose normal direction is tilted from that in the smooth case, 
specified by local zenith and azimuth angles 0 and <p. The slope of a facet along two 
orthogonal directions that are perpendicular to the local zenith direction, say x and y 
directions, can be specified by 

Zx=^ = (H' 2 -l) ,/2 cos(p 

ox , (3a) 

Z y ~ — -(p -1) sintp 

& , (3b) 

where p=cos0. The derivatives in Eqs.(3a) and (3b) are confined to the facet. Since there is 
no sufficient quantitative experimental information on the roughness of ice crystal surfaces 
at present, we assume that irregular roughness of ice crystal surfaces is similar to that of a 
wavy sea surface, which can be specified by the Gram-Charlie distribution [Cox and Munk, 
1954]. If the tilt distribution of the roughness is azimuthally homogeneous (i.e., 
independent of angle <p ), the statistical probability density function for the condition that the 
slopes of a facet along the two axis directions are given by the first order Gram-Charlie, or a 
2-D Gaussian distribution, as follows: 

P(Z x ,Z y ) = ^ expf-CZ^ + Z 2 y )/o 2 ], (4) 

where cr is a parameter determining the magnitude of roughness. It is noted that <7=0-0.005, 
0.005-0.05, 0.05-0.3 correspond to slight, moderate, and deep roughness in the single- 
scattering calculation, respectively. The technical details in accounting for the surface 
roughness in the light scattering computation on the basis of geometric optics method have 
been presented by Yang and Liou [1998]. 
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c. Radiative Transfer Model 


Radiative transfer calculations for cirrus are performed using the adding/doubling 
method. The adding/doubling principle has been expressed mathematically in a matrix form 
[Twomey, 1966; Hunt and Grant, 1966] and in an integral form [Hansen and Travis, 
1974]. A concise formulation in a discrete form for the adding/doubling method is provided 
in Appendix, where we also address some numerical concerns in the radiative transfer 
computation such as the truncation of forward peak of the phase function and a stable 
expansion of the phase function in terms of the re-normalized Legendre function. The 
discrete expression of the adding/doubling principle is straightforward and efficient in 
numerical realization. The present adding/doubling computational program has been 
validated with respect to the various cases presented by Lenoble [1985] and also in 
comparison with DISORT [Stamnes et al., 1988, 2000] for a number of canonical 
problems. 

3. Development of Cirrus Models 

a. The Cirrus 3-Layer Model 

The vertically inhomogenous nature of cirrus clouds was observed during the FIRE- 
II program. Figures la and lb show two different vertical profiles based on replicator 
images of ice crystals in cirrus clouds collected on November 25 and December 5, 1991 
[Heymsfield et al. http://box.mmm.ucar.edu/science/cirrus ]. Three distinct regimes of ice 
crystals are evident from the replicator data. In the uppermost layer, small nonspherical 
“quasi-spheres” are predominant. These particles are essentially nonspherical shapes with 
an aspect ratio approaching 1 . The term “quasi-spherical” used in the analysis of observed 
data and in many theoretical studies is often misleading because the optical properties of 
spherical and nonspherical particles are significantly different even if the nonsphericity of 
particle geometry is not substantial. The middle layer of cirrus is primarily composed of 
pristine ice crystals with well-defined hexagonal shapes or bullet rosettes. The bottom layer 
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contains larger but irregular ice crystals. Roughness can also be noted from the replicator 
images of the irregular ice crystals. In both images, it becomes apparent that the particles 
increase in size and the shapes become more complex from the top to the base of the cirrus. 

The dominant habits of ice crystals and size distributions have been observed with a 
fine vertical-resolution for the two cases shown in Fig.la and lb (28 and 33 layer data were 
obtained in FIRE-II in situ observations for the November 25 and December 5 cases, 
respectively). Based on these data, we constructed the percentages of the various ice crystal 
habits and size distributions for the top, middle, and bottom layers of cirrus clouds. For a 
given layer centered at z with thickness of Az, the mean size distribution and percentage of a 
specific habit are given by 

1 fZ+Az/2 

n{D,z) = —\ n{D,z>)dz>, (5a) 

AzJz-Az! 2 


J >z+Az/2 

z-Az/2 


. _ J z-Az/2 
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Jz-Az/2 


n(D,z’)dz l 


where D is the maximum dimension of ice crystals. The preceding averaging procedure is 
also applied to obtain the mean size distribution and habit percentage for cirrus by assuming 
a single size distribution to represent the cirrus cloud. 

To illustrate the differences between the use of a 3-layer model and a single-layer 
cirrus model, we first confine our discussion to the November 25 case. The thicknesses of 
the top, middle, and bottom layers are 0.79, 0.73, and 1.18 km, respectively. The size 
distributions and crystal habit percentages for the three layers are shown in panels A-C of 
Fig.2a. For the uppermost cirrus layer, “quasi-spherical” habits are dominant for the small 
particle mode (D<100 pm), whereas bullet rosettes are dominant for the large particle mode 
(D>100 pm). The middle layer is composed primarily of columns and bullets. In the bottom 
layer, ice crystals are mainly irregular aggregates, but the edges seem to be rounded, perhaps 
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due to the effects of sublimation. The percentages of various ice crystal habits for the entire 
cirrus layer with a thickness of 2.7 km are shown in the panel D of Fig.2a. It can be noted 
that two geometries (sphere and hexagon) are used for the small ice crystals, as illustrated in 
panels A and D of Fig.2a. The two particle morphologies are used for understanding the 
sensitivity of the radiative properties of cirrus to these small ice crystals, which is of concern 
in remote sensing and radiative transfer modeling. This is because these small ice crystals 
are often misidentified as spheres in observations based on the particle images with blurred 
edges. Even using an optical microscope, the shapes of small ice crystals are unlikely to be 
seen clearly because of the poor instrumental resolving power caused by optical diffraction 
phenomenon [ Ohtake , 1970]. 

Shown in Fig2b are the size distributions and habit percentages for the case of 
December 5. The thicknesses of the three layers (1.24, 1.12, and 1.17 km) in this case are 
more uniform in comparison with the case of November 25. In addition, the former is 
thicker. The December 5 case had a much colder cloud top and the crystal population was 
more pristine, in that many columns and few polycrystals (rosettes) were observed. The 
division of large and small modes for ice crystal habits is at 50 pm for the case of 
December 5. It can be noted that for the small mode, “quasi-spherical” ice crystals are 
dominant. There are substantially large numbers of “quasi-spherical” particles even for the 
middle layer. According to Figs. 2a and 2b, both size distributions and ice crystal habits vary 
substantially with altitude. It should be pointed out that a three-layer or multiple-layer cirrus 
model is realistic as a general rule, as is evident from replicator data shown and also from 
the composite cirrus morphology based on radar observation [Mace et al., 1997]. In the 
following discussions, we select the December 5 case to investigate the effect of vertical 
inhomogeneity on single- and multiple-scattering properties of cirrus. 

b. Radiative properties of cirrus layers 

We employ the scattering computational model described in Section 2 to compute 
the extinction cross sections, single-scattering albedos, and phase functions for ice crystals. 
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Ice crystals are assumed to be orientated randomly in the atmosphere. First, to characterize 
the bulk properties of size distribution, we define the mean maximum dimension for a given 
size distribution as follows: 


<D>= 


f Dn(D)dD 



fA- 

n{D)dD 


( 6 ) 


where £> min and D max are the cutoffs of size distribution at small and large sizes, 
respectively. Foot [1988], Francis et al. [1994], Fu [1996], and Wyser and Yang [1998] 
have found that with this definition for mean particle size, the details of the size distribution 
are not important in specifying the bulk optical properties of cirrus. This feature has also 
been observed in the case of water clouds composed of liquid droplets whose scattering 
properties can be solved using Mie theory [Hansen and Travis, 1974; Hu and Stamnes, 
1993]. Following these studies, we define the effective diameter D e and effective radius r e 
for nonspherical ice crystals with a combination of various habits as follows: 


f “ !V i (D)f i (D)n(D)dD 

D e = 2r e =--^-i , (7) 

" I MD)MDMD)dD 
i 


where fi(D) is the percentage of a specific habit at size D. The summation is carried out for 
all the ice crystal habits. We note that the preceding definition of effective radius reduces to 
that defined by Hansen and Travis [1974] in the case of spherical particles, i.e., 

r e =< r 3 >/ < r 2 > . The mean extinction cross section, single-scattering albedo, and phase 
function are given by 
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Figure 3 shows the bulk microphysical and optical properties for the size 
distributions measured for the December 5 case. The upper row shows the geometric 
configuration of the three layers of cirrus and the mean maximum dimension and effective 
size of ice crystals in these layers. The second and third rows show the extinction 
coefficient and single-scattering albedo. The lower two rows provide the asymmetry 
parameter of the phase functions and the fraction of delta transmission [Takano and Liou, 
1989a] in scattered energy. Note that the delta-transmission is an artifact pertaining to the 
ray-tracing technique, which can be circumvented by using a more accurate physical optics 
approach [. Mishchenko and Macke, 1998]. In the present GOM2 calculation based on a 
simplified algorithm [ Yang and Liou, 1996], we do not account for the spreading of the rays 
associated with delta-transmission for size parameter larger than 100. The use of either a 
geometric optics method or a physical optics approach in dealing with delta transmission in 
the single-scattering calculation for large size parameters does not make a significant 
difference in the radiative transfer computation. The mean maximum dimension of the ice 
crystals is 30 pm, 80 pm, and 132 pm for the top, middle, and lower cirrus layers, 
respectively. The mean effective diameters are 47, 92, and 89 pm for the top, middle, and 
bottom layers, respectively. The bottom layer has a smaller ratio of volume to projected area 
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in comparison with the middle layer because of the porous structures of ice crystal 
aggregates. 

The optical properties of ice crystals are computed at two MODIS bands centered at 
0.65 and 2.1 1 pm. In the computation, the MODIS instrumental response functions at these 
channels have been accounted for. At MODIS 0.65 and 2.11 pm bands, the maximum 
extinction coefficient occurs in the middle layer because the number concentration of ice 
crystals is highest for this layer. At the 2.11 pm band, the single-scattering albedo is higher 
in the top layer than in the lower layers due to the prevalence of small crystals. At both two 
bands, the asymmetry parameter is smaller in the top and bottom layers than in the middle 
layer. In the top layer, the asymmetry parameter is reduced because the particles are smaller 
and the forward peak of the phase function is weaker, while in die bottom layer, the reduced 
asymmetry parameter is caused by the particle roughness. The vertical variability of delta 
transmission is similar to that of the asymmetry parameter. In the top layer, delta- 
transmission is substantially reduced due to the ray-spreading effect associated with small 
size parameters [ Yang and Liou, 1996]. In the bottom layer, the roughness of particle 
surface also reduces the delta-transmission effect. 

Figure 4 shows the phase functions associated with the single-scattering properties 
shown in Fig. 3. At the 0.65 pm band, substantial differences between the phase functions 
for the bottom layer and other layers can be noted at scattering angles near 120°. The phase 
function values for the bottom layer are much larger than those for the top and middle layers 
in the scattering region between 5° and 20°. Evidentiy, the phase function values computed 
by assuming that cirrus are homogeneous are significandy different from the phase 
functions of the three layers in some specific scattering regions. At the 2.1 1 pm band, the 
forward peaks of the phase functions are smaller than at the 0.65 pm band because the size 
parameters are smaller. It can also be noted that the magnitude of the phase function in the 
side scattering and backscattering regions are much lower at die 2.1 1 pm band than at 0.65 
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pm band. In addition, the differences of the phase functions for different layers are 
enhanced at the 2. 1 1 pm band because the particle absorption is strongly sensitive to the 
particle sizes. 

4. Results 

A. Comparison of reflectance feature for two cirrus models 

We use the radiative transfer model described in Appendix to compute the bi- 
directional reflectance to compare the radiative features of the three-layer and single-layer 
cirrus models. To understand the physics in the comparison, one needs to interpret the 
scattering geometry involving the sun and a satellite. For a given solar geometry specified 
by ( 6 s ,(p s ) and view geometry specified by (0 v ,<p v ) , the corresponding scattering angle is 
given by 

0 = cos -1 [- cos6^ cos0 v + sin0 5 sin 0 V cos <j >] , (9) 

where <j>= cp s -(p v is the relative azimuth angle between sun and satellite. Note that 9 S and 
9 V are the inclination angles measured from zenith. Figure 5 illustrates the contours of 
scattering angles versus solar- and view-zenith angles for four cases of azimuth angles. The 
solar zenith and viewing zenith angles range between 0° and 60°. It can be seen that the 
scattering angles for the region of view-solar geometry considered in the present study are 
essentially for side scattering and backscattering directions. The variational pattern of the 
scattering angle versus solar-zenith and viewing-zenith angles depends on the relative 
azimuth angle. A similar contour diagram of the scattering angle versus cos(^) and 
cos(6> v ) has been presented by Mishchenko et al. [1996], 

Figure 6 shows the relative difference between the computed bidirectional 
reflectances of the three-layer and one-layer cirrus models for optically thin cirrus (x = 1) at 
0.65 pm. The relative difference is defined as 
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where R$ and R\ are the bidirectional reflection functions computed using the three-layer 
and one-layer cirrus models, respectively. The maximum difference in this case is about 5%. 
When the optical thickness is small, the photons originating from single scattering events 
dominate the total radiance. The contribution of single-scattering events to the radiance in 
the three-layer case is given by 

r(0 s ,G v ,<t>) = - — Y [©,4^(0,, 0 V ,0)], (11) 

4cos0 t cos0 v ^ 
s v /= 1 

where the summation is carried for all three layers of cirrus. Thus, for thin cirrus the 
bidirectional reflectance function is linearly proportional to the phase function. Referring to 
Figs. 3 and 5, the contours shown in Fig.6 can be explained as follows. For azimuth angles 
of 0° and 60°, the maximum difference is observed near scattering angles of 120°, which 
corresponds to the phase function difference at these scattering angles. For azimuth angles 
of 120° and 180°, the maximum difference for the three-layer and one-layer results are 
mainly noted near scattering angles of 155° and 180°. 

Figure 7 is similar to Fig.6, except that the optical thickness of the cloud is 10. The 
contribution of multiple scattering increases with increasing optical thickness. The 
differences occur at scattering angles between 90° and 120°, between 150° and 160°, and 
between 170°- 180°. From Figs. 6 and 7, the difference between the three-layer and one- 
layer models is within a few percent regardless of optical thickness. Based on these results, 
modeling the cirrus as a single layer would seem to be sufficient at visible wavelengths. 

Figures 8 and 9 are similar to Figs. 6 and 7, except that the calculations are 
performed at the 2. 1 1 pm band. Absorption by ice at the 2. 1 1 pm band is much higher than 
at the 0.65 pm band. It is noted that the reflectance for the three-layer model is always 
larger than the one-layer model at the 2.11 pm band. With the absorption of ice at the 2.11 
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Jim band, the top layer of cirrus dominates the contribution to the cloud reflectance because 
photons that penetrate into the lower layers are largely absorbed. The mean size of the ice 
crystals in the top layer is much smaller than that associated with the general one-layer 
cirrus model. With the increase of optical thickness, the difference between the three-layer 
and one-layer models increases. For an optical thickness of 10, the differences can reach up 
to 50%. Because the 2. 1 1 pm band is used for the retrieval of the mean size and optical 
depth of cirrus cloud, it is expected that the vertical inhomogeneity may be important to 
developing more realistic cirrus retrieval algorithms. In comparing Figs. 8 and 9, it may be 
noted that the difference between the three-layer and one-layer models depends mainly on 
the scattering angle when the cloud is optically thin. However, for optically thick cirrus, the 
difference depends not only on scattering angle but also strongly on the viewing zenith and 
solar zenith angles. This is because the radiance path varies with the solar and view angles. 
For large solar zenith or viewing zenith angles, the ray path is large and the difference 
between the three-layer and one-layer cirrus models, and their associated single-scattering 
properties, becomes more significant. 

B. Sensitivity of cirrus reflectance to shapes of “quasi-spherical” particles 

As discussed in Section 3a regarding the replicator images of ice crystals, the small 
so-called “quasi-spherical” ice crystals are essentially not spherical. Whether their 
morphologies can be treated as spheres in light scattering and radiative transfer calculations 
depends on the wavelength of radiation. At a far infrared wavelength, say, 15 pm or larger, 
the nonsphericity of these particles may not be significant in comparison with wavelength 
because the size parameter becomes small. In this case, spherical and spheroidal geometries 
have been used in light scattering computations [e.g., Sun and Shine, 1994; Takano et al., 
1992]. However, for the 0.65 and 2. 1 1 pm bands, the size parameter for a small particle, say 
on the order of 20 pm, is approximately 100 and 30 at the band centers, respectively. For a 
size parameter on these orders, the nonsphericity effect of a particle on its scattering 
properties can be significant. 
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To investigate the sensitivity of cirrus optical properties to the shapes assumed for 
the small “quasi-spherical” ice crystals, in this study we compare the results associated 
with spherical and hexagonal (with an aspect ratio of unity) assumption for the 
morphologies of these particles. The “quasi-spherical” particles primarily affect the top 
and middle layers and are largely absent in the bottom layer. The effective sizes for the top 
and middle layers are 47.5 pm and 93.0 pm if hexagonal shapes are used for the “quasi- 
spherical” particles whereas the sizes are 52.2 pm and 94.6 pm if perfect spherical 
geometry is assumed for these particles. Evidently, the assumption of shape for the “quasi- 
spherical” ice crystals in cirrus clouds can lead to a change of effective size by as much as 
5 pm. 

Table 1 lists the single-scattering properties of the top and middle layers. Substantial 
differences are noted for delta transmissions at MODIS 0.65 pm band because of the 
absence of parallel faces in spheres. The asymmetry factor of the phase function is 
substantially different at the 2.11 pm band, showing the dependence of the scattering 
properties on the assumption of habits. The differences of the results for the two shapes are 
larger for the top layer than for the middle layer because the former has a larger population 
of “quasi-spherical” ice crystals. 

Figure 10 shows the phase functions associated with the single-scattering properties 
listed in Table 1 . Substantial differences of the phase functions can be seen for the top layer 
at both the 0.65 and 2.1 1 pm bands. The overall feature is that the spherical assumption 
leads to larger phase function values for scattering angles between 10° and 45° and lower 
values at side scattering angles. However, the spherical assumption leads to a pronounced 
scattering maximum between 130° and 140°, which corresponds to the rainbow feature of 
ice spheres. Additionally, the spherical assumption leads to larger phase function values 
near 180° (backscattering) at the 2.1 1 pm band. In the middle layer, the differences caused 
by the assumption of habit for the small particles are reduced due to the small number of the 
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small particles in the size distribution. Evidently, using the most realistic particle shapes for 
the small particles in the top layer of cirrus will be crucial to obtaining reliable single- 
scattering properties of the cloud. 

Shown in Figs. 1 1 and 12 are the difference of bidirectional reflectances computed 
using nonspherical hexagonal and spherical geometries in the three-layer for thin ( T=l) and 
thick (t = 10) cirrus at the 2.1 1 pm band. The differences shown in Figs. 1 1 and 12 are 
defined as 

e(6»„0 v ,^) = lOO[ J R sph «9,,6> v ,0)-i? hex (^,0 v ,0)]//? hex (^,0 v ,^), (12) 

where and i? hex indicate the reflection functions associated with spherical and 

hexagonal shapes, respectively, which are assumed for the small “quasi-spherical” ice 
crystals. The maximum differences shown in Fig. 1 1, which correspond to scattering angles 
between 130° and 140°, are caused by the rainbow feature of ice spheres. It can also be 
noted that the assumption of ice spheres leads to an overestimation of reflectance near 1 80° 
(backscattering). As optical thickness increases, the contrast decreases for the rainbow 
feature. However, the enhanced backscattering for the assumption of spheres is still 
noticeable. For optically thick cirrus, Fig. 12 shows that the assumption of sphere for the 
“quasi-spherical” particles leads to an underestimation of cloud reflection at the 2. 1 1 pm 
band except for scattering angles near 180°. 

Figures 13 and 14 are similar to Figs. 1 1 and 12, except that the formers are for the 
0.65 pm band. Again, we see pronounced differences between the results associated with 
the hexagonal and spherical assumptions for the small ice crystals in the uppermost layer. 
The positive maximum near the backscattering peak noted in Figure 11, however, is not 
observed in the results shown in Figure 13 because the phase function value for the spheres 
is less than that of the hexagons at the 0.65 pm band. Additionally at the 0.65 pm band, a 
distinct rainbow feature can be noted in the case of thin cirrus. For the optically thick cirrus, 
the rainbow is blurred due to multiple scattering events occurring within the clouds. Figures 
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11 through 14 illustrate that the shape effect of the small ice crystals at the top layer of 
cirrus is significant at both visible and near infrared wavelengths. 

5. Conclusions 

In this study, we have defined a three-layer cirrus model in terms of ice crystal habit 
and size distribution based on in-situ replicator data acquired during the NASA FIRE-II 
field observation program. We have described a fundamental scattering model and a 
numerically stable radiative transfer model for the computation of the single-scattering 
properties of various ice crystals and the bidirectional reflection of cirrus clouds. 

We have found that the effect of vertically inhomogeneity within cirrus is not 
significant at MODIS 0.65 pm band, a wavelength for which the absorption of ice is 
negligible. However, in comparison with the one-layer cirrus model, a vertically 
inhomogeneous cirrus cloud produces substantially larger reflectance at the 2.1 1 pm band, a 
wavelength for which absorption by ice is important. The increase in reflectance occurs 
because the mean size of ice crystals in the top layer in the three-layer model is smaller than 
in the case of the one-layer model and the total reflected radiance is dominated by the 
contribution from the top layer. For a given optical thickness, die reflectance increases with 
decreasing particle size. 

We also investigated the sensitivity of reflection of cirrus clouds to the particle 
shapes of the “quasi-spherical” ice crystals that have been often assumed to be spheres. 
Numerical results have demonstrated that the bidirectional reflection function of cirrus is 
very sensitive to the shape of these particles at both visible and near-infrared wavelengths. 
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Appendix : A Discrete Expression of Adding/Doubling Principle 

The adding/doubling method is one of most robust approaches to solve the radiative 
transfer equation for multiple scattering events. The standard mathematical expression of 
this method involves various tedious angular integrals although it can be written 
symbolically in a very simple form. In this section we present a discrete form of the method 
by introducing a direct transmitting function. As an improvement over previous 
mathematical expressions from a practical viewpoint, the discrete adding/doubling equations 
are straightforward and more efficient in numerical implementation. In addition, the discrete 
form of the adding/doubling method is more suitable for addressing some numerical 
concerns, such as the numerical singularity of adding/doubling calculation and the 
performances of various quadrature schemes. 

To economize computational cost and memory requirements, we apply a Fourier 
expansion over the azimuth angle for radiance and bidirectional reflection and transmission 
functions: 

M 

r ,r, ‘(+/j.,(p) = Yj‘' r ’ t{m \+il)cosm(p, (Ala) 

m= 0 
M 

,(p ' ) = )cos m((p-(ff), (Alb) 

m= 0 
M 

,<p ' ) = ^r (m) (At,^' )cosm((p - <p' ). (Ale) 

m= 0 

where superscripts i, r, and t indicate incident, reflected, and transmitted intensities, 
respectively; and p and p' are positive with allowable values in [0,1]. The maximum 
number of Fourier terms ( M) required for the convergent solutions of Eqs.(Ala)-(Alc) 
depends on the incoming and outgoing radiation beams as well as the magnitude of the 
asymmetry of the phase function. For example, only the first term is necessary for a solar 
zenith angle of 0° (overhead Sun) whereas more than 100 terms may be required for a low 
solar elevation angle (Sun near the horizon). For solar zenith and viewing zenith angles 
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smaller than 60°, 30 to 40 terms are normally required in the case of cirrus clouds if the 
strong forward peak of the phase function is truncated. A comprehensive study regarding 
the number of the terms required in the Fourier series has been carried out by King [1983] 
using Henyey-Greenstein function and the phase function for a fair weather cumulus. 

According to the definitions of reflection and transmission functions [. Hansen and 
Travis, 1974], it can be proven that the Fourier components of reflected and transmitted 
intensities are given by 

I r{m) (H) = (l + <5 m0 )f r (m \n,n' )I Km \-it )/i’ dpi , (A2a) 

Jo 

I t{m \v) = (1 + s m0 ) ft" W )/ i(m) (-//' ¥ df £ , (Alb) 

Jo 

where <5 m0 is Kronecker delta function. It should be pointed out that the radiances defined 
in Eqs.(A2a) and (A2b) are diffusive intensities, i.e., they originate from the scattering of 
incident radiation by the particles in the scattering layer. If the total (direct + diffuse) 

transmitted intensity, indicated as I (~p,cp), is defined in the same manner as in the 
diffusive case, it follows that 

I t(m) (H) = (1 + 8 m0 ) f ? (m) (ju,ju' ¥ dfd 

Jo 

= (1 + <5 m0 ) f [t (w W ) + W )]I i(m) (-Li ¥ dii , (A3 a) 

Jo 

where ) is the mth Fourier component of the total transmission function. 

A( m W) is associated with the transmission of incident radiation and given by 

A (m) (/A// ) = - - 1 — exp(-r/ n' )8(n -n'), (A3b) 

(1 + 8 m0 )n’ 

in which 5(p -)i' ) is the Dirac delta function. Evidently, the quantities defined in 
Eqs.(A2a)-(A3b) are continuous functions of the arguments p and p' that range 
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continuously in the interval [0,1]. We select a set of discrete points in the region [0,1] for p 
and p' , that is, the two variables can only have the values of [jU|,p 2 >"Pj • With respect to 

the discrete set, we define the following discrete quantities with single subscript and double 
subscripts: 


if m) = r/ m) = r (m) (fij) , 

r (m) = ( (m) = 


(A4a) 

(A4b) 

(A4c) 


Since the continuous region [0, 1 ] is discretized by using a set of points, it is required that an 
integration of a function f(p) with respect to its argument defined in [0,1] be replaced by a 
discrete summation in the form of 


r 1 ” 

(A5) 
i - 1 

where W t are the weights in the summation. For an integral involving the Dirac delta 
function, the definition of the delta function and Eq.(A5) lead to the following relationship: 


f 1 * ^ 

f(Mi) = - ViW -* ^fQijWjSQij - &)■ (A6) 

j = 1 

Evidently, to guarantee the equality in Eq.(A6) in the discrete procedure, the Dirac delta 
function should be replaced by Kronecker symbol in the form of 


8(n J -H i )^>8 jl /Wj = 



for j = i 
for j*i ' 


(A7) 


Thus, the direct transmission function in discrete form is given by 



1 

(1 +S m0 )W jt ij 


eX P (-T /flj)8y . 


(AS) 
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For the discrete quantities defined with respect to the set of discrete points [ p l5 p 2 > — Mu 1 , we 
introduce a mathematical operator <8> defined by 


n 


4“' ® = 0 + S„ 0 )£ 4" >S S? Vy V , , 

7=1 

(A9a) 

c<" , ®4r > =a+^o 

(A9b) 


7=1 


The operator ® is similar to an ordinary matrix multiplication except that a weight is 
included in the former. Thus, for one homogeneous layer, the reflected, diffusely 
transmitted, and total transmitted radiances are related to the incident radiation via the 
following relationships: 

7 r(m) = jiim) 0 r jm) > jt(m) = ji(m) 0 ,(*), Jt{m) = ji{m) 0 ? jm) ( A 10) 

One of the interesting features of using the operator 0 is the variation of subscripts in the 
expressions in Eq.(A10): the incident beam denoted by subscript j is re-directed to the 
direction denoted by subscript k after interacting with the scattering layer. Similarly, for two 
layers indicated by superscripts a and b, we have the following relationships: 

I r k (m) = lf m > ® [rj? w) + tjf m) ® c4 w) ® ( A1 la ) 

and 

V k {m) = lf m) <g> ® F* (m) + f/ m) 9 ® i$ m \ (A1 lb) 

where the asterisk indicates that the transmission function corresponds to the case of 
illumination coming from below. The quantities D and U in Eqs.(Alla) and (A lib) are 
given by 

N 

D\f = , Ujf = r b(m) + Dj m) 0 r b{m \ (A12a) 

n=l 
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£)(jn),(n+\) _ £){m),(n) ^ £)(m),(l) _ ^(m) ^ r a*(m) . 


(A 12b) 


According to the sensitivity study by Hansen and Travis [1974], we use N=12, 5, and 3 in 
the summation involved in the first expression in Eq.(A13a) for m<10, 10<m<100, and 
m>100, respectively. The remaining terms are approximated by a geometric series. The 
physics of the adding/doubling principle can be viewed clearly in terms of the variations of 
the subscripts from left to right in the right-hand sides of Eqs.(Alla) and (A lib). The 
reflection and transmission functions for the combined layer are given by 

/$"> = r$ m) + ?jf m > 0 Uj^ 0 /£<«>, (A1 3a) 

and 

ij k m) = tjf m) 0 tg m) + 0 0 t$ m \ (A 13b) 

The transmission function given in Eq. (A 13b) contains the contribution due to direct 
transmission, which is implicitly in the form of a delta function, and it may potentially cause 
inaccuracy in numerical computation. Thus, it is necessary to separate the diffusive and 
direct components in Eq.(A13b). It is noted that 

?/"> ® = [<;/“> + 47V,,)] ® [4 ( ”> + 4?>(i t )i = 

® 4 ( "> + exp(-T +<;, ( ">ex p(-l t /ft)+4f> (t„ + T b ). (A14) 

Thus, the diffusive transmission function for the combined layer is given by 
7 <”> = $">»(*<”> 

+lp exp(-T„ //!*) + exp(-r„ //l,)!#"' +?/"> ® ® '5"'- ( A15 ) 

Eqs. (A 12a), (A12b), (A13a) and (A15) constitute the adding/doubling equations. It should 
be pointed that the matrix associated with direct transmission is diagonal. In numerical 
computation, the numerical efficiency and accuracy can be improved substantially if the 
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multiplication of the direct transmission function with another quantity is evaluated 
analytically, such as 

4 m) (T) ® A ( jf = exp(-T/A ® Ag\z) = exp(-r I fi k ). (A16) 

As mentioned in preceding discussions, the continuous integration involved in 
adding/doubling calculation must be replaced by a proper quadrature scheme. 
Mathematically, we need to select proper pairs of (jlpW;) for the definition of the 
mathematical operator ® in Eqs. (A9a) and (A9b). The most popular quadrature schemes 
are Gauss, Lobatto, and Radau quadrature schemes [Press et al., 1986; Hildebrand , 1974]. 
For the angular region involved, these three schemes cover p (or ji') e (0,1), [0,1], and 
(0,1], respectively. That is, Gaussian quadrature is open at both ends, Labotto quadrature is 
closed at both 0 and 1, and the Radau scheme is open at 0 but closed at 1 . The radiance data 
at nadir view is usually required in retrieval applications. Thus, an extrapolation must be 
used to obtain the ndair-view radiance if Gaussian quadrature is used. Although the Labotto 
scheme is closed at both ends of the integral region, the information at p=0 actually does 
not make any contribution to the angular integration, as evident from Eqs.(A9a) and (A9b). 
In addition, including p=0 will cause a singularity in the initialization of the adding/doubling 
calculation. Therefore, we use the Radau scheme in the present study. 

To initialize the adding/doubling process, we start with a very thin homogeneous 

layer with optical depth At (~ 10 -8 ). The reflection and transmission truncation for this 
layer can be obtained on the basis of the invariance principle [Hansen and Travis, 1974] as 
follows: 

(m) _ &AT^ , At Pj + /f i_ pr(m) 

iJ 4 Wj [ 2 HiHj J iJ 

DLfliflj 
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Am) _ 
ij 


d>A T 


4 MiMj 


1 Ax A + ^ 

2 fitfij 


pt(m) 

‘J 


a> 2 A t 2 
32/i,/iy 


[(4 < " , //'i)®(<" > ^*) + (^S < ” ) //</i)®(/S < " ) ^t)]- (A 17b) 


In Eqs.(A17a) and (A 17b), and are defined as follows 


M M 

/J<"> = ^(2-8„ 0 )®,^'(-n,)^”(n,)=^(2-8„„)cb,^"(n,)^"(-n,), (A 18a) 

\~m l=m 

M 

rf"" = £(2 - (A 18b) 

l=m 

in which &>/ is given by 


G), = [ PinWtiLiW 

2 J-\ 

iPM+l-D'PhmivW, (A19) 

2 Jo 

where P(ji) is phase function and P t (p) is the Legendre polynomial of /th order. In 
numerical computation, we use the Radau quadrature scheme for the integration in 
Eq.(A19). Thus phase function information at the exact forward and backward directions is 

accounted for. The function P™ in Eqs.(A18a) and (A18b) are the “re-normalized” ( or the 
“normalized” called in DISORT [Stamnes et al. 2000]) associated Legendre polynomials 
first introduced by Dave and Armstrong [1970], defined as 

<A20a) 

v (i + m ) ! 


where pf is the ordinary Legendre function. The normalized associated Legendre functions 
can be calculated on the basis of the following recurrence relationship: 
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2/ + 1 


~J(l + m + + 1) 


(A21a) 

V (/ + m + 1)(/ -m + 1) 


with the two initial values for the preceding recurrence given by 


tfOO = (-!)” 


(2m-l)! l, 
1 (2m)!! 


OV)”' 2 , 


(A21b) 


Note that alternatives for initializing the preceding recurrence can be found in the paper by 
Dave and Armstrong [1970] and a technical report for DISORT [Stamnes et al. 2000], It 
should be pointed out that in many references the phase function expansion based on 
Legendre polynomials is given in the form of 

M 

= £(2 - (A22a) 

l-m 

where the associated Legendre polynomial or Legendre functions can be calculated on the 
basis of the following recursive relationship 


= (A22b) 

/-m + 1 

The factor (/-m)!/(/ + m)! in Eq.(A22a) rapidly reduces to zero while the values of 
Legendre functions are very large for a large m with / > m, as noted by Dave and 
Armstrong [1970]. For this reason, the preceding approach given by Eqs.(A22a) and 
(A22b) for the expansion of the phase function in radiative transfer simulations is not 
numerically stable, in particular, when the asymmetry of phase function is substantial and 
higher order Legendre functions are required in the phase function expansion. 

Since the predominant sizes of ice crystals in cirrus clouds are much larger than 
visible and near infrared wavelengths, there is a strong forward peak in the corresponding 
phase function. To include this forward peak in numerical computations, thousands of terms 
may be required in the Fourier expansions involved in Eqs.(Ala)-(Alc) for a general solar- 
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view geometry, and also in the phase function expansion given by Eqs. (A 18a) and (A 18b). 
It should be pointed that the number of radiance streams used in the adding/doubling 
calculation needs to increase with the increase of the terms used in the phase function 
expansion so that the orthogonality of the selected finite set of Legendre functions in a 
discrete form can be guaranteed. Thus, the strong forward peak makes the numerical 
computation impractical. For this reason, the strong forward peaks of phase functions are 
truncated in practice. The schemes used for the truncation are diverse, which have been 
evaluated in an extensive discussion by Wiscombe, [1977], who has further developed the 
5 — M method in order to avoid the shortcomings pertaining to various ad hoc empirical 
approaches. For truncating the forward peak of the phase function involved in this study, we 
employ the method developed by Hu et al. [2000], which is an extension of the 5 - M 
method. After the phase function is truncated, die single-scattering properties need to be 
adjusted on the basis of the similarity principle [Wiscombe, 1977; van de Hulst, 1980; 
Takano and Liou, 1989b]. For example, optical depth and single scattering albedo are 
adjusted as follows: 

r = (i-yS>)T, ©’=(i-/)©/(i-y©), (A23) 

where f is the fraction of energy associated with the truncated forward peak. 
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Figure Captions 


Figure la. Replicator images of ice crystals from a cirrus cloud observed on Nov. 25, 1991 
during the FIRE-II field experiment. Note the three-layer structure with small 
quasi-spherical crystals in the top layer, and columns and bullet rosettes in the 
second layer. The third layer is composed mostly of large aggregated crystals. 

Figure lb. Same as Figure la, except that the observation was made on Dec. 5, 1991 and the 
top layer is dominated by pristine columns. 

Figure 2a. Size distribution modeling the cirrus observed on Nov. 25, 1991 that is shown in 
Figure la. Panels A, B, and C show the size and habit distributions for the top, 
middle, and bottom layer, respectively. Panel D shows the mean size distribution 
averaged over height. 

Figure 2b. Same as Figure 2a except for the case of Dec.5, 1991. 

Figure 3. The mean size and single-scattering properties for the three-layer (Dec. 5, FIRE- 
II) cirrus model. The vertical lines indicate the results computed using the one-layer 
model mean size distribution (i.e., the cloud is assumed to be vertically 
homogeneous). 

Figure 4. Phase function corresponding to the single-scattering properties shown in Figure 
3. 

Figure 5. The scattering angle versus solar zenith and view zenith angles for four azimuthal 
angles. Note that the scattering angles are essentially for side-scattering and 
backscattering directions. 

Figure 6. The percent relative difference of bidirectional reflectance computed using the 
three- and one-layer models at MODIS 0.65 pm band for thin cirrus (t= 1). The 
maximum difference for this case is about 5% and depends mainly on scattering 
angle. 

Figure 7. Same as Figure 6 except for thick cirrus (x=10). 
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Figure 8. Same as Figure 6 except for MODIS 2.1 1 pm band. Note the relative difference 
is much higher (up to 12%) than at 0.65 pm wavelength due to absorption by ice. 

Figure 9. Same as Figure 7 except for MODIS 2.11 pm band. The differences for large 
optical thickness reach up to 50% and depend also on viewing and solar zenith 
angles. 

Figure 10. Comparison of the top and middle layer phase functions computed by assuming 
that the small “quasi-spherical” ice crystals are either spheres or non-spherical 
hexagons with an aspect ratio of unity. Note the presence of the ice sphere rainbow 
feature between 130° and 140°. 

Figure 11. The percent relative difference of the bidirectional reflectances computed 
assuming spherical and hexagonal shapes for the small “quasi-spherical” ice 
crystals. The difference contours shown are for thin cirrus (x=l) at MODIS 2.1 1 
pm band. Note the large differences at the ice rainbow and backscattering angles. 

Figure 12. Same as Figure 1 1 except for thick cirrus (t= 1 0). Note the smoothing of the 
rainbow maximum. 

Figure 13. Same as Figure 1 1 except for MODIS 0.65 pm band. Note the absence of the 
positive backscattering angle maximum. 

Figure 14. Same as Figure 13 except for thick cirrus (x=10). 
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Table 1. The single scattering properties of the top and middle layers, which are 
computed in conjunction with two assumptions for the shapes of the “quasi- 
spherical” small ice crystals ( ice spheres and hexagons with aspect ratio of 1). 
Note that the shapes of ice crystals with size larger than 50 pm are assumed to 
be unchanged . 
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